package nl.unimaas.micc.eventSnookerTest.math;

import java.util.ArrayList;
import java.util.Arrays;

public class BairstowSolver {
	
	private double[] coefficients;
	private ArrayList<ComplexNumber> roots = new ArrayList<ComplexNumber>();
	private int order;
	private double threshold = 1e-6;
	
	/**
	 * Polynomial solver using Bairstow's method
	 * @param coefficients 
	 * 				an array with the coefficients stored in the way:
	 * 				f(x) = coeff[0] + coeff[1]*(X^2) + ... + coeff[order]*(X^order)
	 * @param order
	 * 				the order of the polynomial
	 */
	public BairstowSolver(double[] coefficients, int order)
	{
		this.setCoefficients(coefficients);
		this.setOrder(order);
	}

	public void setCoefficients(double[] coefficients) {
		//round too small coefficients so that (hopefully) no infinite loop occurs
		if(coefficients != null){
			for(int i = 0; i < coefficients.length; i++){
				if(Math.abs(coefficients[i]) < 1e-20){
					coefficients[i] = 0;
				}
			}
			
		}
		this.coefficients = coefficients;
	}

	public double[] getCoefficients() {
		return coefficients;
	}

	public void setRoots(ArrayList<ComplexNumber> roots) {
		this.roots = roots;
	}

	public ArrayList<ComplexNumber> getRoots() {
		return roots;
	}

	public void setOrder(int order) {
		this.order = order;
	}

	public int getOrder() {
		return order;
	}
	
	public boolean foundReal()
	{
		for (ComplexNumber root : roots)
		{
			if (root.isRealOnly())
			{
				return true;
			}
		}
		return false;
	}
	
	public double getMinReal()
	{
		if (!foundReal())
		{
			throw new RuntimeException("No real root found");
		}
		double minReal = Double.MAX_VALUE;
		for (ComplexNumber root : roots)
		{
			if (root.isRealOnly() && root.getReal() > 0)
			{
				minReal = Math.min(minReal, root.getReal());
			}
		}
		return minReal;
	}
	
	/**
	 * Solves the given polynomial and puts the results in the arrayList
	 * @return true if there are any roots, false otherwise
	 * @author Stephen R. Schmitt
	 * 
	 *  translated into Java from http://home.att.net/~srschmitt/linbairstow.html
	 *  
	 */
	
	public boolean solve() {
//		System.out.println("BS start");
		roots.clear();
		double[] a,b,d;
		double alpha1 = 0, alpha2 = 0;
		double beta1 = 0, beta2 = 0;
		double delta1 = 0, delta2 = 0, delta3 = 0;
		int i,j,k = 0;
		int count = 0;
		int n;

//		long time1 = System.nanoTime();
		n = this.getOrder();
		
		b = new double[n+2];
		d = new double[n+2];
		

		if (Math.abs(coefficients[0]) < threshold)
			return false;

		a = new double[n+2];
		for(i = 0; i <= n; i++){
			a[n+1 - i] = coefficients[i];
		}
		if (Math.abs(a[1] - 1) > threshold){
			for(i = 3; i <= n+2; i++)
				a[i-1] = a[i-1] / a[1];
			a[1] = 1;

		}
		
		count = 1;
		
		do{
			//guesses
			alpha1 = a[n-1]/a[n];
			beta1 = a[n-2]/a[n];
			do{
				b[0] = 0;
				d[0] = 0;
				b[1] = 1;
				d[1] = 1;
				
				j = 2;
				k = 1;
				
				for(i = 3; i <= n+2; i++){
					b[i-1] = a[i-1] - alpha1 * b[j-1] - beta1 * b[k-1];
					d[i - 1] = b[i-1] -  alpha1 * d[j - 1] - beta1 * d[k-1];
					j = j + 1;
					k = k+1;
				}
				j = n;
				k = n - 1;
				
				delta1 = (Math.pow(d[j-1],2) - (d[n] - b[n]) * d[k-1]);
				if(delta1 == 0)
					return false;
				alpha2 = (b[n] * d[j-1] - b[n+1] * d[k-1]) / delta1;
				beta2 = (b[n+1] * d[j-1] - (d[n] - b[n]) *b[n] ) / delta1;
				
				alpha1 = alpha1 + alpha2;
				beta1 = beta1 + beta2;
//				System.out.println(alpha1);
//				System.out.println(beta1);
//				System.out.println(alpha2);
//				System.out.println(beta2);
//				System.out.println("Threshold: " + threshold);
				
			} while(Math.abs(alpha2) > threshold && Math.abs(beta2) > threshold);
			
			
			delta1 = Math.pow(alpha1,2) - 4 * beta1;
			
			//roots are imaginary
			if(delta1 < 0){
				delta2 = Math.sqrt(Math.abs(delta1)) / 2;
				delta3 = -alpha1 / 2;
				ComplexNumber solution1 = new ComplexNumber(delta3,delta2);
				ComplexNumber solution2 = new ComplexNumber(delta3,-delta2);
				roots.add(solution1);
				roots.add(solution2);
				//System.out.println(solution1);
				//System.out.println(solution2);
			}
			
			//roots are real
			else{
				delta2 = Math.sqrt(delta1);
				ComplexNumber solution1 = new ComplexNumber((delta2 - alpha1)/2,0);
				ComplexNumber solution2 = new ComplexNumber((delta2 + alpha1)/-2,0);
				roots.add(solution1);
				roots.add(solution2);
				//System.out.println(solution1);
				//System.out.println(solution2);
			}
			count = count + 2;
			n = n - 2;
			if(n >= 2){
				for( i = 2; i <= n+2; i++){
					a[i-1] = b[i-1];
				}
			}
			
		} while(n >= 2);
		
		if(n == 1){
			ComplexNumber solution1 = new ComplexNumber(-b[2],0);
			roots.add(solution1);
			//System.out.println(solution1);
		}
	
//		long time2 = System.nanoTime();
//		System.out.println("Running time: " + (time2-time1) + "ns");
//		System.out.println("BS stop");
		return true;
	}
	
	public String toString(){
		String s = "Coefficients:  \n" + Arrays.toString(coefficients) + "\nRoots:\n";
		for(ComplexNumber c : roots)
			s = s + c + "\n";
		return s;
	}
}
